function [moments,residuals] = Var1Threshold_GMM_moments_macro(theta2,VarU_x2,...
    Cov10U_x2,Cov01U_x2,x2,macros,econAct,thresholdLevel,h0RegimeOn,hxRegimeOn)

% Extending the output from th regression filter with macro variables
nm          = size(macros,1);
[nx2,T]     = size(x2);
factors     = [macros;x2];
nx          = size(factors,1);
VarU        = zeros(nx,nx,T);
Cov10U      = zeros(nx,nx,T);
Cov01U      = zeros(nx,nx,T);
for t=1:T
    VarU(nm+1:nm+nx2,nm+1:nm+nx2,t)   = VarU_x2(:,:,t);
    Cov10U(nm+1:nm+nx2,nm+1:nm+nx2,t) = Cov10U_x2(:,:,t);
    Cov01U(nm+1:nm+nx2,nm+1:nm+nx2,t) = Cov01U_x2(:,:,t);
end

% We unfold theta2
[h0_1,h0_2,hx_1,hx_2,sigma,gama0,gamaz,gamax,sigmaz] = unfoldTheta2_threshold_macro(theta2,nx);
VarW          = sigma*sigma';      
regime_1      = econAct >= thresholdLevel; %a boom
regime_2      = econAct < thresholdLevel;  %a recression
     
% The moments conditions 
residuals   = zeros(nx,T);
residuals_1 = zeros(nx,T);
residuals_2 = zeros(nx,T);
residuals_z = zeros(1,T);
tmpx        = zeros(nx,nx,T);
tmpx_1      = zeros(nx,nx,T);
tmpx_2      = zeros(nx,nx,T);
tmpx_w      = zeros(nx,nx,T);
tmpz        = zeros(2+nx+1,T);
for t=1:T-1
   residuals(:,t+1) = factors(:,t+1) - (h0_1 + hx_1*factors(:,t))*regime_1(t,1)...
                                     - (h0_2 + hx_2*factors(:,t))*regime_2(t,1);
   tmpx(:,:,t+1)    = residuals(:,t+1)*factors(:,t)' ...
                        - Cov10U(:,:,t+1) + hx_1*VarU(:,:,t);          
   % For the process for x
   residuals_1(:,t+1) = residuals(:,t+1)*regime_1(t,1);
   residuals_2(:,t+1) = residuals(:,t+1)*regime_2(t,1);
   tmpx_1(:,:,t+1)    = residuals(:,t+1)*regime_1(t,1)*factors(:,t)' ...
                        - Cov10U(:,:,t+1)*regime_1(t,1) + hx_1*VarU(:,:,t)*regime_1(t,1);
   tmpx_2(:,:,t+1)    = residuals(:,t+1)*regime_2(t,1)*factors(:,t)' ...
                        - Cov10U(:,:,t+1)*regime_2(t,1) + hx_2*VarU(:,:,t)*regime_2(t,1);   
   omegaHat = VarU(:,:,t+1)...
       + hx_1*VarU(:,:,t)*regime_1(t,1)*hx_1' ...
       + hx_2*VarU(:,:,t)*regime_2(t,1)*hx_2' ...
       - Cov10U(:,:,t+1)*regime_1(t,1)*hx_1' - hx_1*Cov01U(:,:,t+1)*regime_1(t,1) ...
       - Cov10U(:,:,t+1)*regime_2(t,1)*hx_2' - hx_2*Cov01U(:,:,t+1)*regime_2(t,1);
   tmpx_w(:,:,t+1)    = residuals(:,t+1)*residuals(:,t+1)' - VarW - omegaHat;
   
   residuals_z(1,t+1) = econAct(t+1,1) - gama0 - gamaz*econAct(t,1) - gamax'*factors(:,t);
   % For the process for z
   tmpz(1,t+1)        = residuals_z(1,t+1);
   tmpz(2,t+1)        = residuals_z(1,t+1)*econAct(t,1);
   tmpz(2+1:2+nx,t+1) = reshape(residuals_z(1,t+1)*factors(:,t)'+gamax'*VarU(:,:,t),nx,1);
   tmpz(2+nx+1,t+1)   = residuals_z(1,t+1)^2 - sigmaz^2 - gamax'*VarU(:,:,t)*gamax;
end

% The moment conditions 
namesTheta2 = fieldnames(theta2);
numTheta2 = size(namesTheta2,1);
if hxRegimeOn == 1
    hxRegimeNum = 2*nx^2;
elseif hxRegimeOn == 0    
    hxRegimeNum = nx^2;
end
moments = zeros(numTheta2,T-1);
if h0RegimeOn == 1
    moments(1:nx,:)         = residuals_1(:,2:T);
    idx                     = nx;
    moments(idx+1:idx+nx,:) = residuals_2(:,2:T);
    idx                     = 2*nx;
else
    moments(1:nx,:)         = residuals(:,2:T);
    idx                     = nx;
end
for t=1:T-1
    if hxRegimeOn == 1
        moments(idx+1:idx+nx^2,t)           = reshape(tmpx_1(:,:,t+1),nx*nx,1);
        moments(idx+nx^2+1:idx+nx^2+nx^2,t) = reshape(tmpx_2(:,:,t+1),nx*nx,1);
    else
        moments(idx+1:idx+nx^2,t)        = reshape(tmpx(:,:,t+1),nx*nx,1);
    end
end
idx = idx + hxRegimeNum;
for t=1:T-1
    moments(idx+1:idx+nx*(nx+1)/2,t) = vech(tmpx_w(:,:,t+1));
    moments(idx+nx*(nx+1)/2+1:idx+nx*(nx+1)/2+2+nx+1,t) = tmpz(:,t+1);
end

end